Perfect simulation of spatial point processes 
using dominated coupling from the past with 
application to a multiscale area-interaction 
point process 



We consider perfect simulation algorithms for locally stable point pro- 
cesses based on dominated coupling from the past. A version of the al- 
gorithm is developed which is feasible for processes which are neither 
purely attractive nor purely repulsive. Such processes include multiscale 
area-interaction processes, which are capable of modelling point patterns 
whose clustering structure varies across scales. We prove correctness of 
the algorithm and existence of these processes. An application to the 
redwood seedlings data is discussed. 

1 Introduction 

One of the long standing problems in Markov chain Monte Carlo is that it is 
rarely possible to know when the Markov chain we are using for simulation has 
reached equilibrium. For certain classes of problem, this problem was solved 
by the introduction of coupling from the past (CFTP) [16j[17]. More recently, 
methods based on CFTP have been developed for perfect simulation of spatial 
point process models (see for example [T51 ITUl IT3] ). 

Exact CFTP methods are therefore attractive, as one does not need to rigor- 
ously check convergence or worry about burn-in, or use complicated methods to 
find appropriate standard errors for Monte Carlo estimates based on correlated 
samples. Independent and identically distributed samples are now available, so 
estimation reduces to the simplest case. Unfortunately, this simplicity comes at 
a price. These methods are notorious for taking a long time to return just one 
exact sample and are often difficult to code, leading many to give up and return 
to nonexact methods. 

In response to these issues, in the first part of this paper we present a 
dominated CFTP algorithm for the simulation of locally stable point processes 



Graeme K. Ambler 
University of Bristol 



Bernard W. Silverman 
University of Oxford 



2004 



Abstract 



which potentially requires far fewer evaluations per iteration than the existing 
method in the literature [T3] . 

It is often the case that advancements in statistical theory are inspired by 
applications, and this is no exception. There are several classes of model for 
stochastic point processes, for example simple Poisson processes, cluster pro- 
cesses such as Cox processes, and processes defined as the stationary distribution 
of Markov point processes, such as Strauss processes [H] and area-interaction 
processes [3]. 

All of the above mentioned point process models are capable of modelling 
either clustered or regular point patterns. They are not, however, well suited 
to modelling point patterns whose clustering structure varies across scales, for 
example clusters of regularly spaced points or regularly spaced clusters of points. 
In the second part of the paper we introduce a new multiscale area-interaction 
process which is capable of modelling cither of these types of point pattern. We 
then demonstrate how the algorithm developed in the first part of the paper 
may be used to generate samples from this process. 

The structure of this papers is as follows. In Section [2] we discuss perfect 
simulation, beginning with ordinary coupling from the past (CFTP) and moving 
on to dominated CFTP for spatial point processes. We then introduce and jus- 
tify our perfect simulation algorithm. In Section [3] we first review the standard 
area-interaction process. We then introduce our multiscale process, describe 
how to use our new perfect simulation algorithm to simulate from it, and dis- 
cuss a method for inferring the parameter values from data. An application to 
the Redwood seedlings data is presented in Section [4] and some areas for future 
work are discussed in Section [U 

2 Perfect simulation 
2.1 Coupling from the past 

The principle behind CFTP is the following. Suppose that it is desirable to 
sample from the stationary distribution of an ergodic Markov chain {Z t } on 
some (finite) state space X with states l,...,n. It is clear that if it were 
possible to go back an infinite amount in time, start the chain running (in state 
Z_oo) and then return to the present, the chain would (with probability 1) be 
in its stationary distribution when one returned to the present (i.e. Zq ~ tt, 
where tt is the stationary distribution of the chain) . 

Suppose now that we were to set not one, but n chains {Z^ 1 '}, . . . , {Z^} 
running at a fixed time — M in the past, where Z^ M — i for each chain {Z^}. 
Now let all the chains be coupled so that if Zg = Z& at any time s then 
Z\ ' = Z\ JI Vi > s. Then if all the chains ended up in the same state j at time 
zero (i.e. Z^ = j Vi G X), we would know that whichever state the chain 
passing from time minus infinity to zero was in at time —M, the chain would 
end up in state j at time zero. Thus j must be a sample from the stationary 
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distribution of the Markov chain in question. 

When performing CFTP, a useful property of the coupling chosen is that it 
be stochastically monotone as in the following definition. 

Definition 1 Let {Z t } and {Z^ } be two Markov chains obeying the same 
transition kernel. Then a coupling of these Markov chains is stochastically 
monotone with respect to a partial ordering < if whenever z[ l%> < zj: , then 
zl+ k < ^t+k f or a ^ positive k. 

Whenever the coupling used is stochastically monotone and there are maxi- 
mal and minimal elements with respect to < then we need only simulate chains 
which start in the top and bottom states, since chains starting in all other states 
are sandwiched by these two. This is an important ingredient of the dominated 
coupling from the past algorithm introduced in the next section. 

Although attempts have been made to generalise CFTP to continuous state 



spaces (notably [TS] and [9], as well as [2], discussed in Section 2.2), there 
is still much work to be done before exact sampling becomes universally, or 
even generally applicable. For example, there are no truly general methods for 
processes in high, or even moderate, dimensions. 



2.2 Dominated coupling from the past 

Dominated coupling from the past was introduced as an extension of coupling 
from the past which allowed the simulation of the area-interaction process [15] . 
though it was soon extended to other types of point processes and more general 
spaces [14]. We give the formulation for locally stable point processes. 

Suppose that we wish to obtain a sample of a spatial point process with 
density / with respect to the unit rate Poisson process, whose Papangelou con- 
ditional intensity, \f(u;X), is uniformly bounded above by some constant A: 

/(iuM) 

Xf{u;X) = — jjj^ — < A - 

The uniform bound on the Papangelou conditional intensity is required in order 
for the point process to be locally stable, so this is not imposing any additional 
constraint. 

Then the algorithm given in [14] is as follows. 

1. Obtain a sample of the Poisson process with rate A. 

2. Evolve a Markov process D(T) backwards until some fixed time — T, using 
a birth-and-death process with death rate equal to 1 and birth rate equal 
to A. The configuration generated in step 1 is used as the initial state. 

3. Mark all of the points in the process with U[0,1] marks. We refer to the 
mark of point x as P(x). 
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4. Recursively define upper and lower processes, U and L as follows. The 
initial configurations at time — T for the processes are 

U- T (-T) = {x : x £ D{-T)} 
L_ T (-T) - {0} 

5. Evolve the processes forwards in time to t = in the following way. 

Suppose that the processes have been generated up a given time, u, and 
suppose that the next birth or death to occur after that time happens at 
time ti. If a birth happens next then we accept the birth of the point x 
in U—t or L-t if the point's mark, P(x), is less than 

xmn{^^:L_ T (t i )CXCU_ T (t i )] or 

max^iMtijCIC^)} (1) 

respectively, where x is the point to be born. 

If, however, a death happens next then if the event is present in either of 
our processes we remove the dying event, setting J7_x(ti) = U-t( u ) \ {x} 
and L-x(ti) — L-t{v) \ {x}- 

6. Define U-t{u + e) = U-t(u) and L-t(u + e) = L^t(u) for u < u + e < ti- 

7. If U-t and L-t are identical at time zero (i.e. if U-t(0) = L-t(O)), then 
we have the required sample from the area-interaction process with rate 
parameter A and attraction parameter 7. If not, go to step 2 and repeat, 
extending the underlying Poisson process back to —(T+S) and generating 
additional U[0, 1] marks (keeping the ones already generated). 

This algorithm involves calculation of X(u; X) for each configuration that is 
both a subset of U(T) and a superset of L(T). Since calculation of X(u;X) is 
typically expensive, this calculation may be very costly. The method proposed 
in Section |2.3| uses an alternative version of step [5] which only requires us to 
calculate A(u; X) for upper and lower processes. 

The more general form given in |14j may be obtained from the above algo- 
rithm by replacing the evolving Poisson process D(T) with a general dominating 
process on a partially ordered space (fi, ■<) with a unique minimal element 0. 
The partial ordering in the above algorithm is that induced by the subset rela- 
tion C. Step [5] is replaced by any step which preserves the sandwiching relations 

L^ T (u) r< X^ T {u) r< U- T (u) r< D(u) and (2) 
L_ T (t) = E7_ T (f) if L_ T (s) = U^t(s) (3) 

for s < t < 0, and the funnelling property 

L- T (u) r< L_ [T+S) {u) r< U_ {T+3) {u) < U- T (u) (4) 

for all k < and T,S > 0. In equation ([2]), X_t(w) is the Markov chain or 
process from whose stationary distribution we wish to sample. 
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2.3 A new perfect simulation algorithm 

Suppose that we wish to sample from a locally stable point process with density 

m 

p(X) = al[f i (X), (5) 

i=i 

where a € (0, oo) and fa : — ► M. are positive valued functions which are 
monotonic with respect to the partial ordering ^ induced by the subset relatiorj^J 
and have uniformly bounded Papangelou conditional intensity: 

/,(xuM) 

A/ * (w;x)= m - K - 

Then clearly 

m 

\ p (u; x) < A = Y[ max ^fi { x \ X) (6) 

2—1 

for all u and x, and A is finite. Thus we may use the algorithm in Section [2. 2| to 
simulate from this process using a Poisson process with rate A as the dominating 
process. 

However, as previously mentioned, calculation of X p (u; x) is typically expen- 
sive, increasing at least linearly in n(x). Thus to calculate the expressions in 
display (JIJ, we must in general perform 2™( f/ - T ( ti ))~™( L - T ( ti )) of these calcula- 
tions, making the algorithm non-polynomial. In practice it is clearly not feasible 
to use this algorithm in all but the most trivial of cases, so we must look for 
some way to reduce the computational burden in step [5] of the algorithm. 

This can be done by replacing step [5] with the following alternative. 

[5j Evolve the processes forwards in time to t — in the following way. 

Suppose that the processes have been generated up a given time, it, and 
suppose that the next birth or death to occur after that time happens at 
time tj. If a birth happens next then we accept the birth of the point x 
in U-t or L^t if the point's mark, P(x), is less than 

n^ 1 [max{A /i (u;[/(T)) J A /4 ( U ;i(T))}/A] or (7) 
UT=i [min {X fi (u; U(T)),X fi (u; L(T))} /A] (8) 

respectively, where x is the point to be born. 

If, however, a death happens next then if the event is present in either of 
our processes we remove the dying event, setting U—x{ti) = U_t{ u ) \ { x } 
and L-t(U) = L—t(u) \ {x}. 

Lemma 1 Step [5|' obeys properties ^ and and is thus a valid domi- 
nated coupling from the past algorithm. 

lr That is, configurations x and y satisfy x ^ y if a; C y. 
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Proof Property ^ follows by noting that 

(§ < \ p (u;X) < (|7f < 1. 

Property ^ is trivial. Property Q follows from the monotonicity of the fi's. 
□ □ 

Theorem 1 Suppose that we wish to simulate from a locally stable point process 
whose density p(X) with respect to the unit rate Poisson process is representable 
in form |75p. Then by replacing Step [5] by Step [5]' it is possible to bound the 
necessary number of calculations of A p (it; X) per iteration in the dominated 
coupling from the past algorithm independently ofn(X). 

Proof Step [5] clearly involves only a constant number of calculations, so by 
Lemma [TJ above and Theorem 2.1 of [H] the result holds. □ □ 

In the case where it is possible to write p(X) in form ^ with m = 1, Step[5j is 
identical to Step [5] This is the case for models which are either purely attractive 
or purely repulsive, such as the standard area-interaction process discussed in 
Section [XT] It is not the case for the multiscale process discussed in Section [372] 
or the model studied in [JJ. 

The proof of Theorem 2.1 in [T3] does not require that the initial configu- 
ration of L-t be the minimal clement 0, only that it be constructed in such a 
way as properties ([2]), ([3| and Q are satisfied. Thus we may refine our method 
further by modifying step [4] so that the initial configuration of is given by 



L_ T (-T) = \x G D{-T) : P(x) < JJ 



min Xf^x'jX) /A 

X,{x} 



(9) 



which clearly satisfies the necessary requirements. 



3 Area-interaction processes 

3.1 Standard area-interaction process 

In the standard case, the area-interaction process has density 

p(X) - a \ N ( x h- m{XeG) 



(10) 



with respect to the unit rate Poisson process, where a is a normalising constant, 
A > is the rate parameter, N(X) is the number of points in the configuration 
X, 7 > is the clustering parameter, G is some compact set in M. d and X © G 
is Minkowski addition: 



X ffiG = {a eR d : a = 



where x £ X and z E G}. 



Here < 7 < 1 is the repulsive case, while 7 > 1 is the attractive case. The 
case 7 = 1 reduces the a homogeneous Poisson process with rate A. 
Figure [TJ gives an example of the construction when G is a disc. 
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Figure 1: An example of some events together with circular "grains" G. The 
events in the above diagram would be the actual members of the process. The 
circles around them are to show what the set X © G would look like. If 7 were 
large, the point configuration on the right would be favoured, whereas if 7 were 
small, the configuration on the the left would be favoured. 



3.2 A multiscale area-interaction process 

The area-interaction process is a flexible model which allows for a good range of 
models, from regular through total spatial randomness to clustered. Unfortu- 
nately it does not allow for models whose behaviour changes at different resolu- 
tions, for example repulsion at small distances and attraction at large distances. 
Some real-world examples of places where we see this sort of behaviour are the 
distribution of trees on a hillside, or the distribution of zebra in a patch of sa- 
vannah. A physical example of large scale attraction and small scale repulsion is 
the interaction between the strong nuclear force and the electro-magnetic force 
between two oppositely charged particles. The physical laws governing this be- 
haviour arc different from those governing the behaviour of the area-interaction 
class of models, though they may be sufficiently similar so as to provide a useful 
approximation. 

We propose the following model to capture these types of behaviour. 
Definition 2 The multiscale area-interaction process has density 

p{X) = aA A r (^) 7i -™^©G 1 ) 7 -™(xeG 2)j (u) 



where a, A and N(X), are as in equation (10 1; 71 G [l,oo) and 72 £ (0, 1]; and 



G\ and G2 are balls of radius r% and r 2 respectively. 

The process is clearly Markov of range max{r!,r 2 }. If Gi D G2, we will 
have small scale repulsion and large scale attraction. If G\ C G2, we will have 
small scale attraction and large scale repulsion. 



Theorem 2 The density (11) is both measurable and integrable. 



This is a straightforward extension of the proof of Baddeley and van Lieshout 
[3 for the standard area-interaction process. For completeness, it is given in 
Appendix [A] 
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Y- T (u) m G 



[x@G)\{Y_ T (u)s,G) 



(xfflG)\(y_ T («)fflG) 



Figure 2: Another look at Figure[T]with some shading added to show the process 
of simulation. Dark shading shows Y_t(u) © G where Y_t(u) is the state of 
either U or L immediately before we add the new event and G could be either 
G\ or Gi. Light shading shows the amount added if we accept the new event. 
In the configuration on the left, x © G — (x G) \ (Y-t{u) ffi G), so that the 



attractive term in (12 1 or (13 1 will be very small, whereas the repulsive term 



will be large. In the configuration on the right we are adding very little area to 
(Y_ T (u) © G) by adding the event, so the attractive term will be larger and the 
repulsive term will be smaller. 



3.3 Perfect simulation of the multiscale process 



Perfect simulation of the multiscale process (111 is possible using the method 
introduced in Section 2.3 Since (11) is already written as a product of three 
monotonic functions with uniformly bounded Papangelou conditional intensi- 
ties, we need only substitute into equations (|&f(9]) as follows. 

Substituting into equation ^ , we find that the rate of a suitable dominating 
process is 

A 7 - m(G2) . 

The initial configurations of the upper and lower process U and L are then found 
by simulating this process, thinning with a probability of 

-m(Gi) m(G 2 ) 

7i 7 2 



for L. 

As U and L evolve towards time 0, we accept points x in U with probability 

(12) 

and accept events in L whenever 



-m((affiGi)\C/_r(«)ffiGi) m(G 2 )-m((i$G2)\i-TW«G 2 ) 

7i 72 



Figure [2] gives examples of the construction (x © G) \ Y_t{u) © G. 
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3.4 Parametric inference 

We use maximum pseudo- likelihood [5, 6 ,0 [11] to estimate the parameters A, 
7! and 72 . 

As we saw in Section [3. 3 1 the Papangelou conditional intensity of our process 

is 

X(u- X) = ^~™ l (( n ® Gl )\ JSf ® Gl )^- m (( u ® G2 )\ x ® G2 ) 
Thus the pseudo-likelihood equations for this model are 

El f -m((u®Gi)\X@Gi) -m((u<BG 2 )\X<BG 2 ) > 

r — / 7i 7 2 au > 

(14) 

> — ^ — = / m((u© Gi) \ X © Gi)x 

^ C4 71 J A 

K * e>1 ^ -m((«eGi)\XffiGi) -ra((ueG 2 )\XffiG 2 )^ 

(15) 

and 

v-^ m((x l 0G 2 ) \IfflG 2 ) f ,, s~i \ 

x — tl = / m((u©G 2 )\X©G 2 )x 

72 J A 

x * 6yl , -m((u®Gi)\X®Gi) -m((«0Ga)\XeG 2 ) , 

A7 1 7 2 "I) 

(16) 

where we recall that A is an arbitrary subset of the window in which we observe 
the point process. Clearly the main difficulty is in estimating the integrals on 
the right hand side of equations (14 1 to (16 1. This problem may be tackled 
directly [2] by noting that the integral in the log pseudo-likelihood 

logPL(0;X) = logA(x i; X) - / X(u;X)du 
can be approximated by 



~ m 
/ X(u; X)du — ^ A(itj; X)wj, 

J A 1 



where Uj are points in A and Wj are quadrature weights. Using and extending 
an observation made by [4], [2] note that if the set {uj : j = 1, . . . , m} contains 
all the events {xi : i = 1, ...,n(X)}, then the log pseudo-likelihood may be 
approximated by 

m 

logPL(0; X) ~ YJjii log Xj - Xj)Wj, (17) 
4=1 
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Figure 3: Redwood seedlings data. Left: The data, selected by [TB] from a 
larger data set analysed by [T§]- Right: Plot of the L function for the redwood 
seedlings. There seems to be interaction at 3 different scales: (very) small scale 
repulsion followed by attraction at a moderate scale and then repulsion at larger 
scales. 



where Xj = \(X,Uj), yj — Zj/u)j and 



1 if Uj S {xi : i = 1, . . . , n(X)} 
if Uj $ {a?, : i = 1, . . . , n(X)}. 



For a fixed point pattern X the right hand side of (|17j) is equivalent to the 
log likelihood of independent Poisson variables Y% ~ Poisson(Afc) taken with 



weights Wk, so (17) can therefore be maximised using standard software for 
fitting Generalised Linear Models, such as that in R. 

In order to put the estimation procedure above into practice, we must have 
values for n and r^, the radii of G\ and Gi respectively. Following the lead 
of 0, we suggest fitting the model for a variety of values of these "nuisance 
parameters" which do not fit into the exponential family model, and choosing the 
values which maximise the pseudo-likelihood. It may be wise to plot estimates 
of some standard functions such as K and G in order to narrow the search 
somewhat. 



4 Redwood seedlings data 

We take a brief look at a data set which has been much analysed in the literature, 
the Redwood seedlings data first considered by [19]. We examine a subset of 
the original data chosen by [18] and later analysed by [8] among others. The 
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Figure 4: L and T function plots of the redwood seedlings data. Left: L- 
function plots of the data together with simulations of the multiscale model 
with parameters R x = 0.07, R 2 = 0.013, A = 0.118, 71 = 2000 and 72 = lO" 200 . 
Dotted lines give an envelope of 19 simulations of the model, the solid line is the 
redwood seedlings data and the dashed line is the average of the 19 simulations. 
Right: The same for the T function. 



data are plotted in Figure [3] We wish to model this data using the multiscale 
model we have introduced. From an inspection of the estimated if -function 
(right pane in Figure |3| of the data using Ripley's edge correction scheme [H] 
we estimate values of R\ and R2 as 0.07 and 0.013 respectively, giving repulsion 
at small scales and attraction at moderate scales. It also seems that there is 
some repulsion at slightly larger scales, so it may be possible to use R 2 = 0.2 
and to model the large scale interaction rather than the small scale interaction 
as we have chosen. 

Fitting the remaining parameters by eye again, we chose values A = 0.118, 
71 = 2000 and 72 = 10~ 200 . The remarkably small value of 72 was necessary 
because the value of R 2 was also very small. It is clear from these numbers that 
it would be more natural to define 71 and 72 on a logarithmic scale. Figure 
[4] shows K and T function plots for 19 simulations from this model, providing 
approximate 95% Monte-Carlo confidence envelopes for the values of the func- 
tions. It can be seen that on the basis of these functions, the model appears to 
fit the data reasonably well. 

The plots show several things: Firstly that the model fits reasonably well, but 
that it is possible that we chose a value of R\ which was slightly too large. Per- 
haps R\ — 0.06 would have been better. Secondly, it seems that the large scale 
repulsion may be an important factor which should not be ignored. Thirdly, 
in this case we have gained little new information by plotting the T function 
- the third order behaviour of the data seems to be similar in nature to the 
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second order structure. 



5 Discussion and future work 

We have developed a new method for perfect simulation of locally stable point 
processes. The main advantage of our method is that it allows acceptance 
probabilities to be computed in 0(n) instead of 0(2 n ) steps for models which 
are neither purely attractive nor purely repulsive. Because of the exponential 
dependence on n, the algorithm of [14] is not feasible in these situations. 

We have also developed a multiscale area-interaction process which incor- 
porates both repulsion and attraction and given a method of simulating this 
process exactly. In addition, we have described a method of parametric infer- 
ence for model fitting and given a small application to the redwood seedlings 
data US]]. 

The practical work of the paper has focused on two-scale models but it is 
clear that in practice it is possible to extend the work to multiscale models in 
the general sense. For example, the sample L-function of the redwood seedlings 
might, if the sample size were larger, indicate the appropriateness of a three 
scale model 

p(X) = a A JV W 7 r m(JC ® Cl) 7 a - ,n(x ® Ga) 73 m(Jf ® Ga) . (18) 
The proof given in Appendix [A] can easily be extended to show the existence of 



this process, and (18) is also amenable to perfect simulation using the method 
of Section |2.3| Because of the small size of the redwood seedlings data set a 
model of this complexity is not warranted, but the fitting of such models, and 
even higher order multiscale models in appropriate circumstances, would be an 
interesting topic for future research. 
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A Proof of Theorem [2] 

We prove a slightly more general result. 

Definition 3 The generalised multiscale area-interaction process has density 
P {X) = a x"(x) 7 -MUi(x)) 7 -u 2 (u 2 (x))^ (1Q) 



where a, X and N(X), are as in equation (10); 71 G [1, 00) and 72 € (0,1]; 
v\ and vi are Borel regular measures; Z\ and Zi are myopically continuous 
functions and Uj = [J x . eX Zj{x{). 



Theorem 3 The density (19) is both measurable and integrable. 



Proof If 71 = 1 then ( 19 1 is simply the repulsive case of the area-interaction 



process and the result holds. If 72 = 1 then ( 19 1 is simply the attractive case 
of the area-interaction process and again the result holds. We now consider the 
case 71 > 1 and 72 < 1- 

Let t > and consider V = {X : v^U^X)) < t}. We will show that 

V is open in the weak topology with respect to the function Ui, and thus that 
vi(U\(X)) is weakly upper semicontinuous. Since we know that upper or lower 
semicontinuous functions are measurable it is then a short road to showing that 



(19) is measurable. 

Pick X € V. Then since vi is regular there is an open set G C x containing 
Ui(X) such that vi(G) < t as well. Now clearly X has no events in the set 
H = {x £ X ■ Z x (x) n G c yt (f>), and Z X {H) C J = {K e % : K n G c ^ <p}, 
which is a closed set in the myopic topology. Clearly X (J) = H, and since 
Z\ is a myopically continuous function this shows that H is also closed. It is 
now easy to see that W = {Y e fw : N(Yh) = 0} (where Yjj is the restriction 
of the configuration Y to the set H) is open in the weak topology with respect 
to the function Ui, and since V is the union of a collection of sets of the form 
W then V is also open. 

This shows that X — > vi(Ui(X)) is weakly upper semicontinuous. Thus 
the map X — > exp(— V\(U\(X)) log7i) is weakly lower semicontinuous. Thus 
X — > exp(— v\{U\(X)) log 71) is measurable. 

By a similar argument, the map X — > exp(— v^iU^iX)) log 72) is weakly 
upper semicontinuous, and thus X — ► exp(— 1/2(^/2 PO) log 72) is measurable. 



Since is clearly measurable this means that (19 1 is measurable 
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To see that (19 1 is integrable note that 



< vi(Ui(X)) < v(x) < oo and 
< v 2 {U 2 {X)) < v{ X ) < oo. 



Now the function /(X) = A^^ is integrable, since this is simply the Radon- 
Nikodym derivative of the Poisson process with rate A with respect to the unit 



rate Poisson process. Hence (19) is dominated by an integrable function and is 
therefore integrable. In fact this shows the stronger result that the generalised 
area-interaction process measure is uniformly absolutely continuous with respect 
to the A-rate Poisson process measure and so its Radon-Nikodym derivative is 
uniformly bounded. □ □ 
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